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SECTION  1 


SUMMARY 


A 


3 

The  problem  addressed  in  this  study  is  the  estimation  of 
ship  track  parameters  from  a set  of  observations  on  differences 
in  acoustic  arrival  time  between  bottom  mounted  hydrophones. 

The  studied  model  of  this  problem  consisted  of  the  following 
features : 

ihe  ship  track  is  assumed  to  be  a straight 
line  which  is  described  by  the  closest 
point  of  approach  (CPA)  to  a reference 
hydrophone,  time  of  CPA,  ship's  velocity, 
and  the  angle  between  the  track  and  a 

line  passing  through  the  hydrophones^. 

/ 

The  hydrophones  are  arrayed  in  a straight 
line  at  known  locations. 

l^he  acoustic  path  between  the  ship  and 
each  hydrophone  is  assiimed  to  be  a straight 
line. 

The  error  between  a measured  value  of  the 
time  difference  of  arrival  and  the  true 
value  is  due  to  limitations  in  sensor  and 
processing  capability. 


The  primary  goal  of  the  study  was  to  determine  the  error 
in  the  estimated  track  parameters  due  to  the  measurement  errors. 
In  this  error  or  sensitivity  analysis,  the  true  track  parameters 
are  used  in  appropriate  calculations,  so  the  error  results  are 
indicative  of  the  "best"  accuracy  available  from  the  system. 

An  additional  goal  was  to  determine  which  subsets  of  the  entire 
set  of  observations  have  sufficiently  accurate  results.  (That 
is,  which  subsets  gave  estimate  errors  nearly  as  small  as  those 
obtained  from  all  observations.) 
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This  method  of  analysis  is  based  on  the  generalized  least  I 

squares  error  theory.  This  theory  was  extended  to  treat  the 
non-linear  problem  represented  by  the  estimation  of  the  track 
parameter.  The  extension  provided  not  only  the  formalism  for 
the  error  analysis,  but  a technique  for  estimating  the  track 
parameters  for  a given  data  set.  The  error  analysis  formalism 
was  implemented  by  constructing  appropriate  computer  codes. 

Several  cases  of  ship  tracks,  data  sampling  and  hydrophone 
spacing  were  studied.  The  technique  for  estimating  the  track 
parameters  from  a real  data  set  is  not,  as  yet,  implemented 
in  a computer  code  form.  However,  issues  concerning  that  im- 
plementation were  addressed  and  resolved  to  the  point  of  in- 
dicating what  features  of  the  present  code  must  be  modified. 

In  this  initial  study,  a limited  number  of  track  parameter 
cases  were  examined  to  determine  the  salient  features  of  para- 
meter variation  on  error  estimation.  Numerical  calculations 

were  based  on  a measurement  "error"  (standard  deviation)  of 
-2 

10  sec,  corresponding  to  a measurement  bandwidth  of  100  Hz. 

In  most  cases  "data"  out  to  slant  ranges  of  5 n.mi.  were  used 
in  the  results  presented  here. 

Results  of  two  structured  sets  of  cases  are  noteworthy. 

One  set  involved  varying  the  CPA  from  1 n.mi.  to  4 n.mi.  in 
increments  of  1 n.mi,  keeping  all  other  parameters  fixed. 

Results,  using  nearly  the  same  number  of  observations  in  each 
case,  showed  slight  degradation  of  the  ratio  of  CPA  standard 
deviation  to  true  CPA  (from  nearly  2.5%  to  nearly  6%).  In 
another  set  of  cases,  the  CPA  and  other  parameters  were  held 
fixed  while  the  track-array  angle  was  varied  through  0°,  10°, 

30°,  60°  and  90°.  Results  indicated  severe  impairment  of 
CPA  estimation  capability  for  the  parallel  (0°  case)  track. 

Detailed  studies  of  the  10°  case  showed  that  observations 
corresponding  to  track  locations  on  both  sides  of  array  endfire 
are  required  for  good  CPA  estimation  accuracy.  Furthermore, 
observations  near  endfire  may  be  particularly  important. 


fi 
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phenomenon  show  that  the  CPA  standard  deviation  decreases  with 
the  number  of  samples  (N)  as  N~^  so  long  as  N is  larger  than 
about  10  (at  one-minute  intervals).  For  smaller  sample  numbers, 
the  standard  deviation  is  sensitive  to  the  particulars  (i.e., 

C-  ship  position,  interval  between  samples,  and  the  particular 

case  parameters)  of  the  sample  set.  This  result  suggests  that 
if  the  data  acquisition  system  is  sample  number  limited  (to 
few  numbers  of  samples),  judicious  choice  of  the  sample  set 
can  result  in  better  CPA  estimation  accuracy  than  a random 
sample  set. 

The  remaining  sections  of  this  report  are  organized  as 
follows : 

( 

• Section  2 outlines  the  mathematical 
formalism  of  generalized  least  squares 
theory  as  applied  to  linear  and  non-linear 

< problems . 

• Section  3 presents  the  application  of 
the  extended  non-linear  theory  to  the 
problem  of  estimating  the  ship  track 
parameters. 

• Section  4 presents  the  cases  studied  and 
their  results. 

< • Section  5 concludes  the  report  with  an 

identification  of  the  issues  involved 
with  the  treatment  of  real  data  for  ship 
track  parameter  estimation. 

< 


t 
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SECTION  2 


AN  OUTLINE  OF  THE  MATHEMATICAL  METHOD 
Linear  System 

The  method  may  be  called  a generalized  least  squares  esti- 
mation and  may  best  be  described  with  an  initial  discussion  of 
a linear  problem.^  Consider  a set  of  linear  relationships, 


Yi  = E 
^ J=1 


A.  . X. 
IJ  3 


i = 1, N 


(1) 


or  in  matrix  form, 


Y = AX, 


(2) 


where  Y is  a column  vector  of  length  N of  quantities  to  be 
measured,  A is  a N x M matrix  of  known  parameters,  and  X is 
a column  vector  of  length  M of  quantities  to  be  estimated. 

The  constraint  M<N  indicates  an  overdetermined  system  for  the 
estimations  of  X. 


An 


observation  of  the  Y^  variable,  call  it  y^ , in  general. 


contains  an  error  e^ , for  which  we  can  write. 


e.  = Y.  - y. 
1 1 1 


i = 1,  . . .N 


(3) 


or 


or 


e = Y-y, 


e = AX-y, 


(4) 

(5) 


where  the  last  two  expressions  are  matrix  form  equivalents  of 
the  first.  Consider  the  scalar  function 


S = e^e 


N 


(6) 
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C 


c 


c 


c 


T 

where  e is  the  transpose  of  e and  is  the  row  vector  composed 
of  the  set  {e^}.  Clearly  S is  the  sum  of  the  squared  errors 
and  minimizing  it  with  respect  to  the  X^'s  gives  a set  of 
equations  which  can  be  solved  forl^^,  i.e.,  the  vector  which 
minimizes  S.  Expanding  on  this  discussion,  we  have 


S 


(7) 


which  we  minimize  by  taking 


i = 1.  .M 


(8) 


Solving  the  system  of  equations  in  (8),  we  obtain 


(9) 


Recall  that  A is  known  and  y is  observed,  so  with  this  composite 

of  numbers,  equation  (9)  gives  the  procedure  by  which  the  esti- 

/\  /\  1 /V 

mate  X is  obtained.  For  reasons  not  indicated  here  , X as  given 

in  (9)  is  the  best  linear  unbiased  estimator  (BLUE)  of  the 

parameters  X. 

We  backtrack  for  a moment  to  discuss  an  important  point. 

If  the  variances  of  the  measurements  y are  known,  either  by  a 
series  of  trials  or  error  estimates  based  on  sensor  resolution 
capabilities,  the  function  S in  (6)  does  not  fully  exploit  this 

information.  For  example,  if  the  measurement  y, , carries  a 

2 2 2 
large  variance,  , with  respect  to  ^2’  ^ 

would  be  reasonable  to  weigh  ®ore  than  y^^  in  arriving  at 

an  estimate  'x'.  Consider  then  a covariance  matrix  C on  the 

measurements  y.  The  diagonal  elements  of  C are  the  variances 

of  y^^,  y2...yjj.  The  off-diagonal  terms  represent  correlations 

among  the  y's,  e.g.. 
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C . . = Cov  (y . y . ) 


(yj  “ 

where  ^^denotes  "the  expected  value  or  "average” . 


(10) 


convenience,  let 


W = C"^, 


for  which  we  note  that  if 


(11) 


^ '^i  '^ij  (diagonal), 


(12) 


then 


~ l/t^jL  (also  diagonal) 


Now  we  modify  (6)  to  read 


S = e W e . 


(13) 


(14) 


This  form  "weighs"  in  a consistent  manner  the  variance  infor- 
mation on  the  observations  y.  Note  that  C and  W are  N x N 
matrices.  Note  further  that  for  diagonal  W,  S reduces  to 


kH  ■'A 


(15) 


indicating  a simple  weighing  case.  Carrying  out  the  minimization 
(8)  on  the  form  (14),  we  obtain  the  weighted  BLUE  estimator 


^ T -1 

X = (A^  W A) 


D y , 


(a"^  W y) 


(16) 


(17) 
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where 


D = (A^  W A)"^  (a'^W)  . 


(18) 


Again,  everything  in  D is  known,  so  (16)  is  the  prescription 
for  obtaining  the  estimato'*' "x"  from  the  observations  y.  We 

shall  continue  to  work  with  the  forms  (16)  through  (18). 

The  discussion,  so  far,  is  focused  on  obtaining  the 
estimator  of  a set  of  parameters  X.  Note  that  is , in  fact, 
a random  variable.  In  particular,  'X'  has  a mean  value  and  a 
covariance  matrix  associated  with  it.  By  the  structure  given 
above,  but  not  explicitly  shown  there,  the  mean  of  X is  equal 
to  the  population  mean  and  is  not  of  particular  concern  here 
However,  we  are  interested  in  the  covariance  matrix  since  it 
relates  how  "errors"  (variances)  in  the  observations  are 
propagated  to  "errors"  in  the  estimator.  This  is  the  Law  of 
Covariance  Propagation  which  we  now  demonstrate. 

Consider  the  linear  form 
V = AZ 

where  Z is  the  vector  of  observed  random  variables,  A is  the 
matrix  of  known  coefficients  and  V is  the  vector  of  dependent 
random  variables.  Let  P be  the  covariance  matrix  of  the  ob- 
servations , 


P.J  = COV  (ZjZj) 


assuming  zero  mean  for  convenience, 


C 
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Forming 


Cov  (v.vj) 


= <i: 

k 

1 

II 

*ik  ''ji  <^k^> 

= EE 

^ik  ^jl  \l 

k 1 


hen  defining  Q-^  = Cov  (v^v.),  we  have 

Q = A P a'^  (19) 


L 


L 


c; 


This  is  the  Law  of  Covariance  Propagation  and  relates  the 
"errors”  in  Z to  the  "errors"  in  V. 

Applying  (19)  to  the  form  (17)  with  (18),  we  obtain  the 
covariance  matrix,  Q,  of  the  estimator  to  be 

Q = D C 

= (a"^  WA)~^  (20) 

where  W C = 1,  and  standard  matrix  algebra  is  used. 

Several  important  points  about  (20)  should  be  recognized: 

• Q is  not  dependent  on  a particular  set 
of  observations,  y,  only  on  covariances 
of  the  observations.  The  elements  of 
C may  be  selected  from  a priori  estimates 


C 
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of  basic  sensing  system  resolution 
capabilities.  Thus,  (20)  can  be  used 
to  propagate  these  irreducible  errors 
to  the  errors  in  the  estimator. 


• If  C is  diagonal,  indicating  that  the 
measurements  y.  are  uncorrelated,  so 
that  C..  = tj.  then 


d/ai")  6,. 


This  implies  that 


(21) 


. * k'-'ki ‘kj  (“) 

Thus,  even  if  Oj,  = a = constant,  the  Q 
matrix  is  not  diagonal,  indicating  that 
the  elements  of  the  estimator,  , are 
not  independent . 

• Let  f be  a linear  scalar  function  of 
the  parameters  X,  viz., 

f = E b-  X.  (23) 

where  the  {b^}  are  known. 

Then  the  variance  of  f is  given  by 

<./■<!:  Ebj  Xj> 

which  relates  how  the  covariance  of 
the  estimator  is  propagated  to  the 
variance  of  a scalar  random  variable. 


( 
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Non-Linear  System 


Consider  a non-linear  system  of  the  form 


Yi  = g (a^i, 


^ni’  ^1’ 


i = 1, . . .N, 


(25) 


where  there  are  N quantities  Y.  to  be  observed,  and  M quanti- 
ties to  be  estimated.  There  are  n known  parameters  a,  each 
of  which  may  take  on  a different  value  depending  on  the  obser- 
vation index,  i,  hence,  the  double  subscript.  The  function 
g is  assumed  to  be  the  same  for  each  observation.  The  form 
(25)  is  called  the  condition  equation. 

2 

We  shall  linearize  the  system  in  the  following  way. 

Select  a point  X'  (i.e.,  known  values  for  X)  and  expand  the 
function  g in  a Taylor  series  about  that  point. 


Define 


Y ^ = g (a 


Xi-.X^j)  = g 

( {a}  ; X^' , . . . 

M 

+ T.  ^ 

(X^-X'k)  + e, 

k=l  k 

X' 

1 higher  order  terms. 

li  ni ' 

X’i...X’^) 

(26) 


A. 


= i£ 


ij  3X 


jl  X = X' 


(27) 
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so  when  we  have  a set  of  observations  , we  form 

(X.  - X'.)  ^ , (28) 

which  is  a linear  system  for  the  differences 


Zj=Xj-X'j  (29) 

k 

In  words,  is  the  difference  between  the  observed  value  y^ 
and  the  value  obtained  from  the  condition  equation  evaluated 
at  the  assumed  point  X' ; Z is  the  difference  between  the  point 
X to  be  estimated  and  the  assumed  point  X'.  Thus,  (28)  is  in 
the  form 


V = A Z + e.  (30) 

This  is  the  same  form  as  (5)  and  the  whole  linear  analysis 
discussed  above  follows  through.  Interpretations,  however, 
of  the  derived  equations  differ  somewhat  from  the  purely  linear 
case.  In  particular,  consider  (20).  Since  C is  the  covariance 
matrix  of  the  observations  y.  and  the  Y . are  not  random  vari- 
ables, then  C is  also  the  covariance  matrix  on  V.  The  A in 
equation  (5)  corresponds  to  the  A in  (30),  but  of  course,  the 
later  is  evaluated  at  X' . The  Q matrix  in  (20)  is  the  covariance 
on  the  BLUE  estimator  Z,  but  since  X'  is  not  random,  the  matrix 
is  also  the  covariance  of  the  estimated  point  X. 

Now  (20)  propagates  measurement  errors  to  the  errors  in 
the  estimate  X,  for  a given  assumed  point  X' . If  in  fact  X' 
is  chosen  to  be  the  true  point,  as  is  done  in  our  specific 
problem  discussed  below,  the  resulting  covariance  matrix  Q 
represents  the  best  accuracy  possible. 
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One  last  point  requires  discussion.  If  we  are  interested 
not  only  in  an  error  or  sensitivity  analysis,  but  also  a way 
to  estimate  X,  then  the  above  mathematical  structure  must  be 
used  in  the  following  way; 

• Assume  a point  X'  and  evaluate  the 
partials  of  (27)  at  X' . 

• Construct  or  assume  a covariance 
matrix  on  the  observations. 

• Calculate  Q in  (20). 

• Apply  (17)  to  obtain 

Z'  = D V,  where  (31) 

T 

D = Q A W (and  is  completely  known 
for  given  X' ) . 

• Use  the  result  (31)  in  (29)  to  find 


X = 'z'  + X' 


where  everything  on  the  right-hand 
side  is  known. 

• The  X given  in  (32)  is  used  as  X", 

a second  approximation  to  the  assumed 
point . 

• Repeat  the  above  procedure  until  some 
convergence  criterion  (as  yet  unde- 
fined) is  satisfied. 


(32) 
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APPLICATION  OF  THE  METHODOLOGY  TO  A SPECIFIC  CASE 

We  shall  consider  the  problems  of  estimating  ship  track 
parameters  from  a time  sequence  of  time  delay  measurements  for 
an  ocean  bottom  line  of  hydrophones.  Figure  3-1  illustrates 
the  geometry,  which  is  not  drawn  to  relative  scale.  For  sim- 
plicity, we  consider  3 phones,  one  defining  the  origin  and 
the  other  two  along  the  y-axis  at  y^  and  yg.  The  surface  ship 
track  (assumed  to  be  a straight  line)  is  characterized  by  4 
parameters , 

• r , the  horizontal  CPA  to  the  "0"— 

CO 

phone  (one  at  origin) 


• t , the  time  of  CPA 
o 


• V,  the  speed  of  ship 

• b,  the  angle  between  the  array  of 
hydrophones  and  the  ship  track. 

The  track  parameters  are  to  be  estimated  (phone  positions  are 
assumed  known)  from  measurements  of  the  time  delays,  which  we 
now  define.  Let  T^(t)  be  the  difference  in  arrival  time  be- 
tween the  phone  at  y^^  (1st  phone)  and  the  zeroth  phone,  T2(2) 
is  the  arrival  time  difference  between  the  phone  at  yg  and 
the  zeroth  phone.  Since  the  ship  moves  along  the  track,  the 
delays  are  functions  of  time.  Taking  time  to  be  a discrete  set 
of  points  t^,...tjj,  we  have  measurements  of  the  form 

^ (tf)  -= 

'^2  ^^1^  " ^2 

'^1  (^2^  = ^3 

^2  ^^2^ 


'^2  ^^n' 


'^2N 
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(33) 


where  we  have  also  defined  the  vector  Y. 

We  now  relate  T^(t)  to  the  track  parameters  and  known 
quantities.  Let  the  ocean  depth  be  d and  the  surface  hori- 
zontal range  variables  be  defined  on  Figure  3-1;  they  are 
lower  case  Roman  s3nnbols . The  corresponding  3-D  ranges  are 
defined  as  upper  case  S3nnbols,  e.g., 

(34) 

Here,  r^^  is  the  horizontal  CPA  to  phone  i.  Let  r^(t)  be 
the  instantaneous  horizontal  distance  from  ship  to  phone  i 
at  time  t.  We  see  that 


' '•ci  * 


+ y.  sin  b 

Cl  CO  1 


Q)-yi  cos  b 


(35) 


then 


(t) 


= d^  + r^^  (t) 

= ^ ("’co  ^i 


Finally , 


+ V(t-to)  - y^  cos  b 


R.  (t)  - Rq  (t) 


(36) 


(37) 


= |d^  + (r^o  + sin  b)^ 
+ I^V(t-t^)  - y^  cos  b 


(38) 


where  c is  the  speed  of  sound  (assumed  known  and  constant). 
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V# 


«* 

«» 


r. 


( 


( 


Equation  (38)  is  the  condition  equation  we  set  out  to 
obtain.  It  is  in  the  form  of  (25)  where  the  jaj  correspond 
to  d,  y^,  Y2>  ^ 

and 


(39) 


We  need  the  partials  indicated  in  (27).  These  are: 


T 

ill 

3r. 


= ("'co  b)/RJt)  - r„^/R^(t) 


CO 


T 

3 i 
3t_ 


= V"(t-t^)/R^(t)  - 


CO'  o 


T 

3 i 
3V 


V I^V(t-tQ)  - y^  cos  tj/R^(t) 


= (t-t^)  V(t-t^)  - y^^  cos  b 


/Ri(t) 


3^1 

3b 


-V(t-t^)"/R^(t) 


= ( ^co  t»  + V (t-t^)  sin  b 


)/Ri 


(t) 


In  order  to  set  up  the  problem  in  the  form  (20), where  Y 
is  given  by  (33)  and  X by  (39),  we  must  define  the  components 
of  A as 


T 

^ (t  ) 
3X.  '■  n-* 
<3 


i = 1,2,  j = 1,  . . .4 
n = 1 , . . .N 


(41) 


)• 

( 


(40) 
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where  k goes  from  1 to  2N  and  is  given  by 


k = 2n  - i 


Thus,  A is  a 2N  (rows)  by  4 (columns)  matrix. 


(42) 


Finally, 


Y.  = V A,  . X. 

k Z-r  kj  j 


k = 1,  . . .2N 


(43) 


We  shall  base  our  computations  on  the  true  track  parameters, 
which  we  specify  for  each  case  we  examine.  That  is,  the  partials 
in  (41)  are  evaluated  at  the  true  values  of  the  parameters  in 
(39).  The  time  set  {t^}  is  selected  in  a way  to  be  described 
later,  so  the  components  of  A are  known.  We  can  then  use  (20) 
for  the  error  propagation  analysis  where  we  specify  C,  the  co- 
variance  matrix  on  the  observations  of  time  delay.  In  the  cal- 
culations to  be  described,  we  have  chosen 


i,  j=l, . . .2N 


(44) 
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where  a =10  sec  , thus,  all  measurements  are  uncorrelated 
with  respect  to  each  other  and  the  standard  deviation  is  given, 
say,  by  the  time  resolution  in  a 100  Hz  bandwidth. 

Let  P = Q”^  for  convenience.  Then  from  (22)  we  have 


_ 1 
a 


A A 
*ki  *kj 


i .J  = 1,4 


(45) 


In  the  computer  code  a matrix 
1 2N 

®ij  ’ L \i  ' 


(46) 
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is  computed. 


Thus 


Q 


2 

a 


(47) 


from  which  we  note  that  the  covariance  on  the  BLUE  of  the  track 
parameters  scale  in  proportion  to  the  variance  on  the  time  de- 
lay observations. 


Equation  (20)  is  the  focus  of  what  has  been  done  to  date. 

One  other  set  of  calculations  is  also  made.  Using  (24),  we 
propagate  the  covariances  in  the  track  parameters  to  the  variances 
of  the  following  quantities; 


• Ship's  X coordinate  at  each  time  step, 

X (t)  = r^Q  cos  b + V (t-t^)  sin  b 

• Ship's  y coordinate, 

y(t)  = -r sin  b + v (t-t  ) cos  b 

CO  o 

• Aspect  angle  from  bow  to  zeroth  phone 


0(t) 


n 


+ 


( 


• Depression  angle 

it)(t)  = tan~^  |d/r^(t) 

• Two  transmission  loss  expressions, 
single  path; 

TL^(t)  = 20  log  R(t)  + 66 
single  path  from  a surface  dipole; 


(48) 


(49) 


(50) 

(51) 


(52) 


where 


2 


i 

I 


R(t)  = 


+ r 


CO 


V(t-t^) 


(54) 


is  the  slant  range  from  the  ship  to 
the  zeroth  phone. 

In  order  to  use  (24)  to  calculate  the  variances  of  these 
parameters,  each  of  the  above  expressions  must  be  linearized 
by  the  Taylor  series'  expansion  about  the  true  (input)  track 
parameters.  The  resultant  linearization  gives  the  coefficients 
{b^}  in  (24). 

The  computer  code,  developed  to  carry  out  the  required 
calculations  of  the  covariance  matrix  in  (48)  and  other  para- 
meters (e.g..  Equations  48  through  54),  is  discussed  in  the 
Appendix . 

Some  of  the  features  of  the  present  formulation  include; 

(1)  We  have  implemented  (46),  which  assumes 
diagonal  C.  The  track  parameters, 

b,  V,  and  t , are  specified  as  input 
parameters. 

(2)  The  time  points  t^  are  selected  by  first 
choosing  a maximum  ship-to-zeroth-sensor 
range  which  is  feasible  under  physical 
conditions.  Using  the  input  value  of 

r , b,  and  V,  the  total  track  length 

CO 

of  interest  is  determined,  which  is 
then  divided  by  V-T^j^^,,  where  is 

the  selected  time  increment  between  the 
observations.  The  time  points  with 
respect  to  t^  are  then  calculated. 

(3)  On  the  first  run  of  a new  set  of  track 
parameters,  a file  is  created  which 
stores  the  time  steps  and  partial 
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derivatives  (evaluated  at  those  time 
steps) , The  covariance  matrix  Q is 
1 evaluated  for  all  the  time  steps. 

(4)  Upon  subsequent  runs,  the  user  can 
call  the  file,  just  described,  and 
specify,  upon  input,  which  time  steps 
are  to  be  used  in  the  calculation. 
Thus,  the  procedure  in  (2)  is  not  a 
limitation  but  simply  a convenience. 
That  is,  the  user  can  look  at  any 
subset  of  times  (to  approximate  a 
physical  condition)  for  evaluating 
the  appropriate  covariance  matrix 
on  the  track  parameter  estimator. 

It  is  noted  that  with  2 sensors  a 
minimum  of  2 time  steps  is  required 
for  a well  posed  solution  (there  are 
4 track  parameters). 


j 
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SECTION  4 
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CASES  AND  PARTIAL  RESULTS 


Cases 

The  track  parameters  examined  so  far  are  summarized  in 
Table  4-1.  Figure  3-1  identifies  the  geometrical  parameters. 
This  limited  set  of  cases  does  not  represent  a complete  vari- 
ation of  track  parameters.  In  the  first  look,  we  wanted  only- 
some  representative  cases.  Note,  however,  two  structured 
sets ; 

(i)  Cases  4,  10,  11,  9 represent  a cross- 
array track  for  which  the  surface  CPA 
varies  from  1 to  4 n.mi.  All  other 
track  parameters  are  equal. 

(ii)  Cases  5,  8,  6,  7,  4 represent  a set 
for  which  the  CPA  is  constant  at 
1 n.mi.,  but  the  track  angle  with 
respect  to  the  array  line  varies. 

One  primary  interest  in  these  initial  data  was  to  deter- 
mine how  the  variances  of  the  estimator,  for  a given  case, 
change  with  respect  to  the  number  of  time  steps  and  the  par- 
ticular time  steps.  Cases  1,  4,  and  8 were  so  analyzed. 

All  calculated  quantities  have  used  a measurement  variance 

-4  2 

of  10  sec  , consistent  with  an  error  of  .01  sec  corresponding 
to  a measurement  bandwidth  of  100  Hz. 
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Table  4-1.  Track  Parameters  for  Cases 
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In  this  case,  the  zeroth  phone  is  between  the  other  two  and  the  ship  passes 
directly  over  the  zeroth  phone. 


Some  Results 


The  results  of  the  cross  track  set  with  varying  CPA  item 
(i),  are  summarized  in  Table  4-2(a).  Results  of  item  (ii)  are 
summarized  in  Table  4-2 (b).  Note,  for  case  5,  where  the  ship 
traverses  parallel  to  the  array,  the  standard  deviation  on  the 
CPA  estimator  is  larger  than  the  known  CPA.  The  indicated 
dramatic  reduction  in  standard  deviation  with  increasing  angle 
suggests  that  a significant  reduction  in  track  error  can  be 
obtained  by  including  some  measurements  with  ship  near  array 
"endfire." 

This  feature  is  supported  by  a study  of  Case  8.  Here, 
several  sets  of  time  steps  were  chosen  and  the  standard 
deviation  of  the  CPA  estimate  was  plotted  versus  the  number  of 
samples  in  a set.  Figure  4-1  siunmarizes  the  results.  The 
crosses  correspond  to  sets  where  the  points  near  end-fire  are 
eliminated.  The  dots  represent- sets  which  extend  the  same  range 
(which  crosses  end-fire)  but  with  fewer  number  of  time  steps  in 
that  range.  Based  on  an  assumed  law 


the  crosses  have  a = 1.69  and  the  dots  a = .50.  Figure  4-1 
may  be  interpreted  in  the  following  way:  Start  at  the  top  of 
the  solid  line  where  the  number  of  samples  is  high,  giving  a 
low  standard  deviation.  If  one  decreases  the  sample  size, 
moving  downward  along  the  solid  line,  by  removing  samples  in 
a way  which  preserves  those  samples  at  or  near  endfire,  the 
solid  line  is  traversed.  If,  however,  the  sample  removal 
systematically  eliminates  the  steps  at  endfire,  the  transition 
to  the  dashed  line  occurs.  Thus,  for  a given  standard  devia- 
tion on  the  CPA  estimate,  say  .5  N.mi.,  sets  of  points  which 
discriminate  against  endfire  or  crossarray  time  steps  require 
larger  numbers  of  time  steps  than  those  sets  that  do  not 
discriminate . 
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Figure  4-2  illustrates  a property  which  may  exist  for 
• most  tracks  and  data  sets  which  include  sufficiently  large 

numbers  of  observations  for  the  ship  on  both  sides  of  the 
array  line.  On  this  figure,  data  points  were  removed  while 
I keeping  the  track  end  points  fixed.  The  relationship  between 

the  CPA  standard  deviation  and  the  sample  size 


I 

i appears  to  have  ot  = .5,  which  is  expected  from  standard  sta- 

tistical theory  for  uncorrelated  samples.  However,  for  small 
1 i sample  sizes  (say,  under  10),  the  standard  deviation  appears 

[.  I to  be  strongly  dependent  on  the  particulars  (i.e.,  ship  loca- 

i f tion  of  the  sample,  sample  interval,  and  specific  case  param- 

[ I eters)  of  the  sample  set. 


Table  4-2.  Partial  Results  For  Two  Parametric  Variations 


NUMBER  OF  OBSERVATI 


Standard  Deviation  of  CPA  Estimate 
(N.Mi. ) 


Numerical  Results  for  Case  8 


The  present  computer  code  represents  a capability  to 
study  the  effect  of  observational  errors  on  track  estimator 
errors.  The  estimate  of  the  track  parameters  (CPA,  time  of 
CPA,  ship  velocity,  and  the  array-track  angle)  is  not  con- 
structed from  the  time  delay  measurements.  In  the  following 
discussion,  we  outline  some  issues  in  implementing  such  a 
capability,  especially  as  an  extension  of  the  present  code 
configuration.  In  Section  2,  an  iterative  methodology  was 
outlined  (discussion  near  Equation  (31)  ff),  which  can  be 
adapted  for  the  purpose  of  track  parameter  estimation.  Be- 
sides the  implementation  of  this  methodology  in  code,  there 
are  two  other  factors  which  must  be  addressed.  These  are 
data  preparation  and  convergence  criteria  selection. 

Data  Preparation 

Data  may  be  supplied  in  a "raw"  form.  Not  knowing  the 
exact  form  of  much  data  at  this  point,  we  assume  that  a 
processor  can  be  constructed  to  supply  a set  of  time  delay 
observations  between  each  of  M phones  and  a reference  on 
zeroth  phone  at  each  of  N times.  The  times  need  not  be 
equally  spaced. 

Some  processing  of  these  time  delays  may  be  desired 
to  calculate  the  covariance  matrix  C of  these  observations. 

If  such  calculations  turn  out  to  be  infeasible,  C may  be 
estimated  from  bandwidth  considerations.  If  a non-diagonal 
ij  C is  calculated  or  assvimed,  the  present  code  must  be  modified 

(easily)  to  accommodate  the  calculation  of  equation  (20) 
rather  than  the  present  use  of  Equation  (22)  with  = constant. 

The  data  processor  is  a iterative  scheme  requiring  a 
first  guess  for  the  track  parameters.  An  adequate  first 


) 

f 

> 
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guess  may  be  constructed  by  examination  of  the  display  from 
one  beam  when  the  ship  crosses  endfire  and  the  beam  center- 
line.  The  times  of  these  events  plus  range  estimates  from 
broadband  correlation  will  uniquely  determine  an  initial  set 
of  track  parameters.  For  ships  not  crossing  endfire  in 
direct -path  range,  two  beam-center  crossings  are  required. 

This  technique  has  not  yet  been  applied,  and  its  convergence 
criteria  remain  to  be  determined. 

Iteration  and  Convergence 

The  iteration  scheme  is  described  earlier  (cf.  equations 
(31),  ff.).  The  required  calculations  in  equation  (31)  are 
presently  carried  out  in  TRACK  and  subroutine  COM  . These 
codes  contain  I/O  functions  which  should  be  omitted  when  the 
iteration  is  carried  out.  The  iterative  step  (32)  involving 
the  construction  of  the  next  guess  to  the  track  parameters 
may  be  easily  implemented.  The  issue  is,  however,  how  to 
condense  the  essential  calculations  in  TRACK  and  COM  to 
efficiently  cycle  through  the  successive  guesses.  The  present 
code  configuration  may  be  modified  with  sufficient  "flags"  to 
accomplish,  this  cycling.  The  resultant  add-ons  to  the  present 
configuration  may  be  unduly  complicated.  It  is  suggested  that 
the  essential  calculations  of  the  present  code  may  be  incor- 
porated into  a new  structure  which  is  more  "tuned"  to  the 
cycling  function  and  the  method  of  deciding  convergence. 

Convergence  may  be  determined  in  the  following  way: 


» Let  X. 


Time  of  CPA 


Xg  ■ Velocity 
X.  ■ Angle 


These  routines  are  described  in  Appendix  A. 


L 
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• Let  be  the  value  of  X.  after  the 

*3  «3- 

iteration . 


• Compose  the  index 


I = 
n 


(n-1) 


where  Wj^  is  an  assigned  weighting  factor. 

For  example,  if  CPA  accuracy  is  more  im- 
portant than  that  of  other  track  parameters, 
one  may  choose  = 4,  = Wg  = = 1. 

• Compare  I with  a predetermined  "tolerance" 

e,  that  is,  the  first  n for  which  I < e turns 

n 

off  the  iteration. 

Modifications  to  this  criterion  are  possible.  One  can, 
for  example,  use  only  the  CPA  parameter  and  check  if 


X<“>  1 x'”-!' 
v(n-l) 


It  is  suggested  that  the  required  accuracies  for  the  track 
parameters  be  based  on  the  physical  problem. 
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APPENDIX 


AN  OUTLINE  OF  THE  PRESENT  COMPUTER  CODE  CAPABILITY 


Two  separate  programs  (TRACK  and  ?IIN)  with  their  associated 
subroutines  and  functions,  constitute  the  present  computing 
capability  implemented  on  a CDC  6600  machine.  Table  A-1 
illustrates  the  structure  of  TRACK  and  its  subprograms. 

This  program  is  used  for  a first  calculation  of  a new  TRACK 
parameter  case.  Basic  results  are  put  on  a file  and  saved 
for  subsequent  use  by  program  MIN  and  its  subprograms. 

The  following  comments  apply  to  calculations  carried 
out  in  TRACK : 

• Calculations  are  performed  in  Hrs-NMi 
units  until  printout  stage  where  some 
conversions  are  performed.  Thus,  in- 
put values  of  parameters  are  converted 
to  these  units. 

• Immediately,  the  horizontal  CPA  variable 
name  is  changed  to  RT  and  RCO  is  used, 
thereafter  for  the  slant  CPA  given  by 

RCO  = SQRT(D*D  + RT*RT'> 

• The  number  of  time  steps  is  calculated 
from 

. 2*SQRT(RMAX*RMAX  - RC0*RC0) 

V*TINC 

where  now  RCO  is  the  slant  CPA  and  TINC 
is  in  hours. 
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• A data  file  "TAPE  7"  is  made  containing 


pertinent  input  and  the  results  of  the 
derivatives  evaluated  at  each  time  step. 

TAPE  7 must  be  categorized  by  the  appro- 
priate Scope  demand. 

The  input  parameters  for  TRACK  are  indicated  in  Table 

A-2. 
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Table  A-2.  TRACK  Input  Parameters 


Variable 

Units 

Algebraic 

Symbol 

Computer 

Symbol 

Horizontal  CPA 

N.MI 

^co 

RCO 

Ship  Velocity 

KNOTS 

V 

V 

Time  of  CPA 

HRS. 

"o 

TO 

Angle  between 
array  and  track 
lines 

DEGREES 

b 

B 

Ocean  depth 

N.MI 

d 

D 

Number  of  sensors 
(excluding  the 
implied  one  at 
the  origin) 

NSS 

Locations  along 

Y axis  of  the 
sensors 

FEET 

^i 

Y(I) 

Maximum  slant 
range  of  interest 

N.MI 

«m 

RtlAX 

Time  increment 
between  samples 

MIN 

^inc 

TINC 

Flag  for  printing 
time  step  infor- 
mation (—1  printing 
is  suppressed) 

IFL 

Flag  (No  longer 
functional  - 
leave  blank) 

— 

— 

JFL 
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The  structure  of  program  MIN  is  outlined  in  Table  A-3. 
The  following  comments  apply  to  the  running  of  MIN; 

• Only  one  of  the  above  options  can  be 
performed  during  a single  execution. 

• Tape  7,  the  data  file  made  by  TRACK, 
must  be  attached.  The  header  infor- 
mation on  TAPE  7 identifying  the  case 
under  consideration  is  printed  so  that 
handy  reference  is  contained  on  the 
output  of  MIN  runs. 

MIN  requires  input  in  either  of  two  forms.  These  two 
input  options  are: 

• The  user  can  select  the  set  of  time 
steps  to  be  used  in  the  covariance 
matrix  calculations.  The  user  supplies 
the  nvimber  of  time  steps,  JFL,  then  the 

indices  of  the  time  steps,  NS(1), NS 

(JFL).  The  program  can  cycle  through 
any  number  of  sets  by  setting  NJFL.  A 
sample  input  deck  is  illustrated  in 
Table  A-4 . 

• The  user  can  select  a sequence  of  time 
steps  by  identifying  IBEG,  the  time 
step  index  of  the  beginning  of  the 
desired  sequence,  lEND,  the  time  step 
of  the  end,  and  INC,  the  index  incre- 
ment. This  option  is  activated  by  set- 
ting a flag  INTNL  f 0.  The  actual 
value  of  INTNL  causes  program  cycling 
using  multiple  IBEG,  lEND,  INC  values. 

An  illustration  is  shown  in  Table  A-5. 


Table  A-3.  Structure  of  MIN 


I 

I 

i 

I 

I 


I 
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MIN 

• Read  file  of  stored  information  on  case  at  hand. 

• Read  cards  for  processing  parameters. 

• If  a flag  is  set  for  special  time  step  selection, 
call  SORT. 

• Call  COM,  which  performed  the  same  functions  as 
in  the  TRACK  program  package,  but  uses  only  the 
time  step  information  provided  to  MIN. 

• COM  calls  VAR,  which  performs  as  in  the  TRACK 
package,  again  using  only  the  provided  time  step 
information . 


Table  A-4.  Option  1 for  MIN 
CARD  SAMPLE 


NO. 

FORMAT 

VARIABLES 

CONTENTS 

COMMENTS 

1 

215 

NJFL.INTNL 

3 (Blank) 

Three  sets  of  time 
steps  will  be  used. 
INTNL  must  be  blank 
or  zero . 

2 

15 

JFL 

4 

Four  time  steps  are 
used  in  the  first 
pass. 

3 

1615 

NS(I); 

1=1, . . .4 

1,14,16,20 

Time  steps  whose 
indices  are  1,14, 
16,20  are  used  in 
the  calculation  . 

4 

15 

JFL 

6 

Six  time  steps  are 
used  in  the  second 
pass . 

5 

1615 

NS(I) 

1=1, . . .6 

1,2, 3, 4, 5,6 

Indices  of  time  steps 

6 

15 

JFL 

2 

Third  pa.ss  has 

2 time  steps. 

7 

1615 

NS(1)NS(L) 

1,  40 

Time  steps  1 and 

40  are  used. 

Table  A-5,  Option  2 for  MIN 


CARD 

NO. 

FORMAT 

VARIABLES 

SAMPLE 

CONTENTS 

COM!JENTS 

1 

215 

NJFL.INTNL 

(BLANK) ,3 

Any  value  for  NJFL  is 
ignored.  Three 
passes  will  be  per- 
formed. 1 

2 

315 

IBEG.IEND, 

INC 

1,6.1 

All  time  steps  bet- 
ween indices  1 and 

6 (inclusive)  are 
used.  (This  is 
another  way  of  doing 
the  second  pass  in 

Option  1 example.) 

3 

315 

IBEG,  lEND, 
INC 

1.6.2 

Time  steps  1,3,5  are 
used.  Note  trun- 
cation . 

4 

315 

IBEG , lEND , 
INC 

1,18,5 

Time  steps  1,6,11, 

16  are  used. 
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